
l <- readRDS("temp/geocoded_shootings_from_gva.rds") |> 
  filter(date >= "2016-09-08",
         date <= "2017-01-08",
         state == "California") |> 
  filter(id %in% (readRDS("temp/geocoded_shootings_from_gva.rds") |> 
                    filter(mass_shooting,
                           home + dv == 0,
                           vics > 3))$id ## drop if in the home or domestic violence
  )

ballot_initiative <- fread("raw_data/state_g16_sov_data_by_g16_srprec/state_g16_sov_data_by_g16_srprec.csv") |> 
  mutate(share_yes = PR_63_Y / (PR_63_Y + PR_63_N),
         mj_con = PR_64_Y / (PR_64_Y + PR_64_N),
         dp_con = PR_62_Y / (PR_62_Y + PR_62_N),
         pb_con = PR_67_Y / (PR_67_Y + PR_67_N)) |> 
  select(SRPREC_KEY, share_yes, mj_con, TOTVOTE, COUNTY,
         pb_con, dp_con, FIPS, ADDIST, TOTREG)

r12 <- fread("raw_data/state_g12_sov_data_by_g12_srprec/state_g12_sov_data_by_g12_srprec.csv") |> 
  select(SRPREC_KEY, PRSDEM01, PRSREP01,
         PR_36_Y, PR_36_N,
         PR_34_Y, PR_34_N)

############################
bs12 <- fread("raw_data/state_g12_sr_blk_map/state_g12_sr_blk_map.csv") 

r12 <- inner_join(r12, select(bs12, SRPREC_KEY, BLOCK_KEY, PCTSRPREC)) |> 
  select(-SRPREC_KEY) |> 
  mutate(across(c(PRSDEM01, PRSREP01,
                  PR_36_Y, PR_36_N,
                  PR_34_Y, PR_34_N), ~ . * PCTSRPREC / 100)) |> 
  select(BLOCK_KEY, PRSDEM01, PRSREP01,
         PR_36_Y, PR_36_N,
         PR_34_Y, PR_34_N) |> 
  group_by(BLOCK_KEY) |> 
  summarize(across(c(PRSDEM01, PRSREP01,
                     PR_36_Y, PR_36_N,
                     PR_34_Y, PR_34_N), ~ sum(., na.rm = T)))

bs16 <- fread("raw_data/state_g16_sr_blk_map/state_g16_sr_blk_map.csv")

r12 <- inner_join(r12, bs16) |> 
  mutate(across(c(PRSDEM01, PRSREP01,
                  PR_36_Y, PR_36_N,
                  PR_34_Y, PR_34_N), ~ . * PCTBLK / 100)) |> 
  group_by(SRPREC_KEY) |> 
  summarize(share_dem = sum(PRSDEM01, na.rm = T) / sum(PRSDEM01 + PRSREP01, na.rm = T),
            share_yes_pr_36 = sum(PR_36_Y, na.rm = T) / sum(PR_36_Y + PR_36_N, na.rm = T),
            share_yes_pr_34 = sum(PR_34_Y, na.rm = T) / sum(PR_34_Y + PR_34_N, na.rm = T))

ballot_initiative <- left_join(ballot_initiative, r12)

shps <- read_sf("raw_data/srprec_state_g16_v01_shp", "srprec_state_g16_v01")

########################################
sf::sf_use_s2(FALSE)
pings  <- st_as_sf(l, coords = c("lon", "lat"), crs = st_crs(shps))
joined <- st_join(pings, shps[, c("SRPREC_KEY")])
l$SRPREC_KEY <- joined$SRPREC_KEY

########################################

centroids <- data.table(lon = st_coordinates(st_centroid(shps))[, "X"],
                        lat = st_coordinates(st_centroid(shps))[, "Y"])

centroids$SRPREC_KEY <- shps$SRPREC_KEY

out <- rbindlist(lapply(unique(l$id), function(shooting){
  
  
  x <- pointDistance(select(centroids, lon, lat),
                     filter(l, id == shooting) |> 
                       ungroup() |> 
                       select(lon, lat), lonlat = T) * 0.000621371
  
  return(cbind(centroids, x) |> 
           mutate(id = shooting) |> 
           rename(distance = x))
  
}))


out <- left_join(out, select(l, id, SRPREC_KEY) |> 
                   mutate(d2 = 0)) |> 
  mutate(distance = ifelse(is.na(d2), distance, 0)) |> 
  select(-d2) |> 
  # filter(distance < 1) |>
  left_join(select(l, id, date)) |> 
  mutate(pre = date <= "2016-11-08") |> 
  left_join(ballot_initiative) |> 
  group_by(SRPREC_KEY, pre) |> 
  filter(distance == min(distance),
         distance <= 20)

##########################


# ca_file <- fread("M:/Democracy & Justice/democracy/voter_file_data/national/post_2016/VoterMapping--CA--HEADERS--2017-03-25/VoterMapping--CA--HEADERS--03-25-2017-HEADERS.tab",
#                  select = c("General_2016-11-08",
#                             "Residence_Addresses_Latitude",
#                             "Residence_Addresses_Longitude",
#                             "Parties_Description"))
# 
# ca_file <- filter(ca_file, !is.na(Residence_Addresses_Latitude))
# 
# 
# pings  <- SpatialPoints(ca_file[,c('Residence_Addresses_Longitude',
#                                    'Residence_Addresses_Latitude')], proj4string = shps@proj4string)
# 
# 
# ca_file$SRPREC_KEY <- over(pings, shps[shps$SRPREC_KEY %in% out$SRPREC_KEY,])$SRPREC_KEY
# 
# 
# ca_file <- filter(ca_file, !is.na(SRPREC_KEY))
# 
# pct_partisanship <- ca_file |>
#   group_by(SRPREC_KEY) |>
#   summarize(dem = mean(Parties_Description == "Democratic"))
# 
# 
# saveRDS(pct_partisanship, "temp/ca_2016_pct_partisanship.rds")


################################



################################
################################
################################

precinct_demos <- readRDS("temp/pct_demos.rds") |> 
  filter(year == 2016)

shp <- read_sf("raw_data/vest/vest_2016/2016/ca_2016") |> 
  mutate(precinct = paste0("CA", row_number())) |> 
  select(precinct, SRPREC_KEY) |> 
  st_drop_geometry()

precinct_demos <- right_join(precinct_demos, shp)

precinct_demos <- left_join(precinct_demos, readRDS("temp/ca_2016_pct_partisanship.rds"))

out <- left_join(out, precinct_demos)

run_thresh <- function(threshold){
  out <- out |> 
    filter(distance <= threshold) |> 
    group_by(SRPREC_KEY) |> 
    filter(pre == max(pre),
           !is.na(FIPS))
  
  out <- out |> 
    mutate(share_other = 1 - (nh_white + nh_black + latino + asian))
  
  m1 <- feols(share_yes ~ pre + median_income + nh_black + latino + asian + nh_white + 
                some_college + median_age + pop_dens + share_dem + dem + share_yes_pr_36 | FIPS, out,
              cluster = ~ c(id, FIPS))
  
  summary(m1)
  
  m3 <- feols(mj_con ~ pre + median_income + nh_white + nh_black + latino + asian +
                some_college + median_age + pop_dens + share_dem + dem + share_yes_pr_36 | FIPS, out,
              cluster = ~ c(id, FIPS))
  
  summary(m3)
  
  m4 <- feols(dp_con ~ pre + median_income + nh_white + nh_black + latino + asian +
                some_college + median_age + pop_dens + share_dem + dem + share_yes_pr_36 | FIPS, out,
              cluster = ~ c(id, FIPS))
  
  summary(m4)
  
  m5 <- feols(pb_con ~ pre + median_income + nh_white + nh_black + latino + asian +
                some_college + median_age + pop_dens + share_dem + dem + share_yes_pr_36 | FIPS, out,
              cluster = ~ c(id, FIPS))
  
  summary(m5)
  
  
  if(threshold == 0.5){
    
    rows2 <- tribble(~term,          ~m1,  ~m3, ~m4, ~m5,
                     "County Fixed Effects", "Y", "Y", "Y", "Y")
    
    attr(rows2, 'position') <- c(25:25)
    
    
    modelsummary(list("Check" = m1,
                      "Penalty" = m4,
                      "Legalization" = m3,
                      "Plastic Bags" = m5),
                 statistic = "std.error",
                 stars = c("*" = 0.05, "**" = 0.01),
                 coef_map = c("preTRUE" = "Near Shooting Before Election",
                              "nh_black" = "Non-Hispanic Black",
                              "latino" = "Latino",
                              "asian" = "Asian",
                              "nh_white" = "White",
                              "median_income" = "Median Income",
                              "some_college" = "Some College",
                              "median_age" = "Median Age",
                              "pop_dens" = "Population Density",
                              "dem" = "Share of Registered Voters who are Democrats",
                              "share_dem" = "Obama Vote Share, 2012",
                              "share_yes_pr_36" = "Share Supporting Three-Strikes Reform (Prop 36), 2012"),
                 gof_omit = 'DF|Deviance|AIC|BIC|Within|Pseudo|Log|Std|FE',
                 title = "\\label{tab:ca-reg} Support for Ballot Initiatives, CA 2016, Precincts Within 0.5 Miles of Mass Shooting",
                 output = "latex",
                 add_rows = rows2,
                 notes = c("Robust standard errors clustered by mass shooting ID and county."),
                 escape = FALSE) |>
      kable_styling(latex_options = c("scale_down")) |> 
      add_header_above(header = c(" ", "Background", "Abolish Death", "Marijuana", "Ban"),
                       line = F) |> 
      save_kable(paste0(outpaths_tabs[2], "ca_reg_full.tex"))
  }
  
  plot_d <- rbindlist(lapply(list(m1, m3, m4, m5), function(m){
    cbind(confint(m), lab = rownames(confint(m)))
  })) |> 
    mutate(p = (`2.5 %` + `97.5 %`)/2,
           t = threshold)
  
  plot_d$model <- c(rep("Prop 63\n(Gun Control)", length(m1$coefficients)),
                    rep("Prop 62\n(Repeal Death Penalty)", length(m4$coefficients)),
                    rep("Prop 64\n(Recreational Marijuana)", length(m3$coefficients)),
                    rep("Prop 67\n(Outlaw Plastic Bags)", length(m5$coefficients)))
  plot_d$model <- factor(plot_d$model, levels = (c("Prop 63\n(Gun Control)",
                                                   "Prop 62\n(Repeal Death Penalty)",
                                                   "Prop 64\n(Recreational Marijuana)",
                                                   "Prop 67\n(Outlaw Plastic Bags)")))
  plot_d$nobs <- m1$nobs
  
  return(plot_d)
}


plot_data <- rbindlist(lapply(c(0.5, seq(1, 10, 1)), run_thresh))

plot_data <- mutate(plot_data, nobs2 = nobs / 80000)

ggplot(filter(plot_data, lab == "preTRUE"),
       aes(x = t, y = p, ymin = `2.5 %`, ymax = `97.5 %`)) +
  geom_point() +
  scale_y_continuous(
    limits = c(-.06, .09),
    # name = "Estiamted Effect Size",                        # Left y-axis label
    
  ) +
  geom_errorbar(width = 0) +
  facet_grid(. ~ model) +
  geom_hline(yintercept = 0) +
  theme_bc() +
  scale_x_continuous(breaks = seq(2, 10, 2)) +
  labs(x = "Radius around shooting (miles)",
       y = "Coefficient and confidence\ninterval on \"Near Shooting Before Election\"")


ggsave(paste0(outpaths[2], "ca_spatial.png"), width = 7.5, height = 5, units = "in")

fwrite(filter(plot_data, lab == "preTRUE") |> 
         select(initiative = model,
                distance = t,
                estimate = p,
                lower_bound = `2.5 %`,
                upper_bound = `97.5 %`), "data_for_figures/figure_A11.csv")
